package org.ainlolcat.ainscience.statistics.handlers.fitting;

import org.ainlolcat.ainscience.statistics.handlers.fitting.functions.ZeroBasedGaussFunction;
import org.ainlolcat.ainscience.statistics.handlers.fitting.functions.ZeroBasedPsdVoigtFunction;
import org.ainlolcat.ainscience.utils.ArraysUtil;
import org.apache.commons.math.FunctionEvaluationException;
import org.apache.commons.math.optimization.OptimizationException;
import org.apache.commons.math.optimization.fitting.CurveFitter;
import org.apache.commons.math.optimization.fitting.ParametricRealFunction;
import org.apache.commons.math.optimization.general.LevenbergMarquardtOptimizer;
import org.apache.log4j.Logger;

/**
 * @author ainlolcat
 *         Date: 10/12/13
 */
public abstract class AbstractPeakFitting implements PeakFitting {

    private static Logger log = Logger.getLogger(ZeroBasedPsvVoigtPeakFitting.class);

    protected ParametricRealFunction function;
    protected double mean;
    protected double disp;
    protected double maximum;
    protected PeakFittingStatisticHandler peakFittingStatisticHandler;

    protected AbstractPeakFitting(PeakFittingStatisticHandler peakFittingStatisticHandler) {
        this.peakFittingStatisticHandler = peakFittingStatisticHandler;
    }

    /**
     * prepare graph for approximation. For example subtracts base line
     * @param x
     * @param y
     * @return
     */
    protected abstract double[] prepareGraph(double[] x, double[] y);

    /**
     * get initial fitting parameters based on mean value (expected as peak center) and dispersion
     * @param mean
     * @param disp
     * @return
     * @throws FunctionEvaluationException
     */
    protected double[] initialGuess(double mean, double disp) throws FunctionEvaluationException {
        double[] initialguess1 = {
                1,                              //A
                0.5,                            //m (m -> L, 1-m -> G)
                Math.abs(disp),                 //w_L
                Math.abs(disp),                 //w_G
                mean                            //x0
        };
        initialguess1[0] = 1/function.value(mean,initialguess1);
        return initialguess1;
    }

    /**
     * is preapproximation supported. for example gauss before psv Voigt
     * @return
     */
    protected boolean isPreApproxSupported(){
        return true;
    }

    /**
     * make approximation parameters prediction
     * @param mean
     * @param disp
     * @param fitter
     * @return
     * @throws OptimizationException
     * @throws FunctionEvaluationException
     */
    protected double[] preapproximation(double mean, double disp, CurveFitter fitter) throws OptimizationException, FunctionEvaluationException {
        double[] initialguess1 = new double[4];
        double[] initialguessG = {1, Math.sqrt(Math.abs(disp)), mean};
        double[] preresult = fitter.fit(function, initialguessG);
        initialguess1[0] = preresult[0];
        initialguess1[1] = 0; //m = > all to Gauss
        initialguess1[2] = 0; //zero width of lorentz
        initialguess1[3] = preresult[1];
        return initialguess1;
    }

    /**
     * get expected width for this approximation
     * @param result
     * @return
     */
    protected double getWidth(double[] result){
        double w_L = result[2];
        double w_G = result[3];
        return (Math.abs(w_L) + Math.abs(w_G));
    }

    /**
     * get center of approximation
     * @param result
     * @return
     */
    protected double getCenter(double[] result){
        return result[result.length-1];
    }

    public double[] getParamsFor(double[] x, double[] y, String name) {
        //prepare graph for fitting. normalize. search for peak center and other
        double[] y1 = prepareGraph(x,y);

        double[][] temp = new double[3][y1.length];
        temp[1] = y1;

        double[][] temp2 = new double[3][y1.length];
        temp2[1] = y1;

        LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer();
        //fill data points
        CurveFitter fitter = new CurveFitter(optimizer);
        for (int i=0;i<(x.length<y.length?x.length:y.length);i++){
            fitter.addObservedPoint(x[i], y1[i]);
        }


        //lets try
        try {
            //prepare initial guess
            double[] initialguess1 = initialGuess(mean, disp);

            //if gauss preapproximation
            if (peakFittingStatisticHandler.gaussPrediction && isPreApproxSupported()){
                initialguess1 = preapproximation(mean, disp, fitter);
            }
            double[] result = fitter.fit(function, initialguess1);

            double x0 = getCenter(result);
            //lets find width
            double max = function.value(x0,result); //value at peak
            double peakWidth = 0;
            double step = getWidth(result)/y1.length;
            for (int i=0;i<200;i++){
                if (function.value(x0+i*step,result)<max/2){
                    peakWidth = 2*(i*step);
                    break;
                }
            }
            for (int i=0;i<y1.length;i++){
                temp[0][i] = x[i];
                temp[2][i]= function.value(x[i], result);
            }
            peakFittingStatisticHandler.data.put(name!=null?name: Double.toString(result[result.length - 1]),temp);
            return new double[] {x0, peakWidth, ArraysUtil.getNormalizedDifference(temp[0], temp[1], temp[2])[1]};
        } catch (Exception e) {
            log.error(e);
            return null;
        }
    }

}
